Protocol 2 - Petrographic data
The following example applies protocol 2 to confirm workshops’ petrographic groups.
Protocol 2 consist in:
- Select ordinal petrographic data;
- Transform to ranks;
- Extended Gower distance, using:
- Relative ranking difference (RRD)
- Neighbor interchange (NI)
- Apply ordination procedure:
- Principal Coordinates Analysis (PCoA)
- Non-metric Dimensional Scaling (NMDS)
- Perform PERMANOVA & PERMDISP tests;
Last, search for outliers and re-do protocol excluding outliers.
NOTE: The initial procedures must be ran at least once before any protocol can be applied.
The key references on the Extended Gower distance are:
Pavoine, S., Vallet, J., Dufour, A.-B., Gachet, S., Daniel, H., 2009. On the challenge of treating various types of variables: application for improving the measurement of functional diversity. Oikos 118, 391-402. doi:10.1111/j.1600-0706.2008.16668.x
Podani, J., 1999. Extending Gower’s General Coefficient of Similarity to Ordinal Characters on JSTOR. Taxon 48, 331-340. doi:10.2307/1224438
Gower, J.C., 1971. A General Coefficient of Similarity and Some of Its Properties. Biometrics 27, 857-871. doi:10.2307/2528823
Ordination procedure
Depending on which type of distance calculation (RRD/NI), protocol 2 performs different ordination methods (PCoA/NMDS). Both PCoA and NMDS require specifying the number of dimensions in which to project the data. Therefore, you must generate specific 2D and 3D ordination objects:
prot2a_2d <- apply_ordination(cleanAmphorae[!isShipwreck,],
"2a", # select protocol 2a (RRD & PCoA)
exception_columns = excep_cols,
variable_tags = varCode)
prot2b_2d <- apply_ordination(cleanAmphorae[!isShipwreck,],
"2b", # select protocol 2a (NI & NMDS)
exception_columns = excep_cols,
variable_tags = varCode)
prot2a_3d <- apply_ordination(cleanAmphorae[!isShipwreck,],
"2a", # select protocol 2a (RRD & PCoA)
exception_columns = excep_cols,
variable_tags = varCode,
dimensions = 3)
prot2b_3d <- apply_ordination(cleanAmphorae[!isShipwreck,],
"2b", # select protocol 2a (NI & NMDS)
exception_columns = excep_cols,
variable_tags = varCode,
dimensions = 3)The ordination objects generated with protocol 2 are different from those in protocol 1 since it uses different functions. However, the main components are still the same: the projection of observations or scores (points) and of variables or loadings.
class(prot2a_2d)
#> [1] "list"
names(prot2a_2d)
#> [1] "points" "eig" "x" "ac"
#> [5] "GOF" "sub2D" "GOF2_2D" "loadings"
#> [9] "variable_tags" "name" "dist_matrix"
class(prot2b_2d)
#> [1] "metaMDS" "monoMDS"
names(prot2b_2d)
#> [1] "nobj" "nfix" "ndim" "ndis"
#> [5] "ngrp" "diss" "iidx" "jidx"
#> [9] "xinit" "istart" "isform" "ities"
#> [13] "iregn" "iscal" "maxits" "sratmx"
#> [17] "strmin" "sfgrmn" "dist" "dhat"
#> [21] "points" "stress" "grstress" "iters"
#> [25] "icause" "call" "model" "distmethod"
#> [29] "distcall" "distance" "converged" "tries"
#> [33] "engine" "species" "data" "init_seed"
#> [37] "trymax" "sub_stress" "sub2D" "GOF2_2D"
#> [41] "loadings" "variable_tags" "name" "dist_matrix"Test the given fabric groups
The fabric groups defined in previous studies can be tested against Protocol 2 distance matrices. We can use either “prot2a_2d$dist_matrix” or “prot2a_3d$dist_matrix”, because they are the same. Remember that this test batch may take several minutes.
prot2a_tests <- test_groups(prot2a_2d$dist_matrix,
factor_list$FabricGroup)
#> [1] "initiating test batch..."
#> [1] "vegan::anosim done."
#> [1] "vegan::betadisper done."
#> [1] "vegan::permutest done."
#> [1] "vegan::adonis done."
#> [1] "Test batch completed."
prot2b_tests <- test_groups(prot2b_2d$dist_matrix,
factor_list$FabricGroup)
#> [1] "initiating test batch..."
#> [1] "vegan::anosim done."
#> [1] "vegan::betadisper done."
#> [1] "vegan::permutest done."
#> [1] "vegan::adonis done."
#> [1] "Test batch completed."These tests were explained in protocol 1.
Biplots
The details on how to create biplots is already explained in protocol 1. Concerning protocol 2, we can compare the results of version 2a (RRD, PCoA) and 2b (NI, NMDS).
Biplot 2D
arrows_label_adj <- rbind(c(.5,.8),c(.5,1),c(.5,1),c(.5,0),c(.5,1),
c(.5,0),c(0,.5))
row.names(arrows_label_adj) <- c("L48","L24","L5","L36","S7",
"S8","S11")
biplot2d3d::biplot_2d(prot2a_2d,
ordination_method = "PCoA",
invert_coordinates = c (TRUE,TRUE),
xlim = c(-.26,.35),
ylim = c(-.31,.35),
point_type = "point",
groups = factor_list$FabricGroup,
group_color = color_list$FabricGroup,
group_label_cex = 0.6,
arrow_mim_dist = 0.5,
arrow_label_cex = 0.6,
arrow_fig = c(.6,.95,0,.35),
arrow_label_adj_override = arrows_label_adj,
subtitle = prot2a_2d$sub2D,
test_text = prot2a_tests$text(prot2a_tests),
test_cex = 0.8,
test_fig = c(0, 0.5, 0.62, .99),
fitAnalysis_fig = c(0,.7,.05,.5),
output_type = "preview")
protocol 2a
arrows_label_adj <- rbind(c(.5,1),c(.5,0),c(.5,1),c(.5,1),c(.5,0),
c(0,.5),c(1,.5))
row.names(arrows_label_adj) <- c("S7","S8","CLAY","L24","L43",
"L5","L36")
biplot2d3d::biplot_2d(prot2b_2d,
ordination_method = "NMDS",
xlim = c(-.42,.38),
ylim = c(-.45,.25),
point_type = "point",
groups = factor_list$FabricGroup,
group_color = color_list$FabricGroup,
group_label_cex = 0.6,
arrow_mim_dist = .5,
arrow_label_cex = 0.6,
arrow_fig = c(.6,.95,0,.35),
arrow_label_adj_override = arrows_label_adj,
subtitle = prot2b_2d$sub2D,
test_text = prot2b_tests$text(prot2b_tests),
test_cex = 0.8,
test_fig = c(0, 0.5, 0.62, .99),
fitAnalysis_stress_axis_cex = 0.8,
fitAnalysis_fig = c(.1,.6,.1,.4),
output_type = "preview")
protocol 2b
Biplot 3D
biplot2d3d::biplot_3d(prot2a_3d,
ordination_method = "PCoA",
point_type = "point",
groups = factor_list$FabricGroup,
group_color = color_list$FabricGroup,
group_representation = "stars",
star_centroid_radius = 0,
star_label_cex = .8,
arrow_min_dist = .5,
arrow_body_length = .025,
subtitle = prot2a_3d$sub3D,
test_text = prot2a_tests$text(prot2a_tests),
test_cex = 1.25,
test_fig = c(0, 0.5, 0.65, .99),
view_zoom = 0.9)
biplot2d3d::animation(directory = directories$prot2,
file_name = "Prot2a_Biplot3D")
Prot2a_Biplot3D.gif
biplot2d3d::biplot_3d(prot2b_3d,
ordination_method = "NMDS",
point_type = "point",
groups = factor_list$FabricGroup,
group_color = color_list$FabricGroup,
group_representation = "stars",
star_centroid_radius = 0,
star_label_cex = .8,
arrow_min_dist = .5,
arrow_body_length = .025,
subtitle = prot2b_3d$sub3D,
test_text = prot2b_tests$text(prot2b_tests),
test_cex = 1.25,
test_fig = c(0, 0.5, 0.65, .99),
view_zoom = 0.9)
biplot2d3d::animation(directory = directories$prot2,
file_name = "Prot2b_Biplot3D")
Prot2b_Biplot3D_snapshot.png